clear;clear;clear;clear;clear;clear;
clc;
%% create session, connect devices
session=MeasurementSession.SimpleMeasurement([Constant.AutoMeasurementRootPath '/example/DemodPhase/device.dat']);

probe_eater_volt = 0:0.025:0.9;
probe_power = [0.594, 0.594, 0.592, 0.592, ...
               0.590, 0.681, 0.825, 0.973, ...
               1.122, 1.271, 1.420, 1.570, ...
               1.719, 1.868, 2.018, 2.167, ...
               2.317, 2.467, 2.617, 2.766, ...
               2.916, 3.065, 3.215, 3.364, ...
               3.514, 3.664, 3.814, 3.963, ...
               4.112, 4.261, 4.411, 4.560, ...
               4.710, 4.859, 5.008, 5.158, ...
               5.307]; % manually measured by power meter
           
result = cell(1, length(probe_eater_volt));
for k = 1:length(probe_eater_volt)
    ziDAQ('setDouble', '/dev3543/auxouts/0/offset', probe_eater_volt(k) );

    %% get Spectrum
    fprintf('\n\nSpectrum: probe power = %f\n', probe_eater_volt(k));
    ch = '3'; lpOrder = 3; lpBW = 100; freqSpan = 400.0; nFreq=2048; average = 100;
    spectrum = getSpectrum('dev3543', ch, lpOrder, lpBW, freqSpan, nFreq, average);
    
    %% store
    result{k}.probe_eaater_volt = probe_eater_volt(k);
    result{k}.spectrum = spectrum;
end

%% plot
for k=1:length(probe_eater_volt)
    plotSpectrum(result{k}.spectrum.dev3543.demods(4).sample_x_fft_abs_avg{1});
    hold on;
end

figure;

bgRange = 100:1000; 
bgNoise=zeros(1, length(probe_eater_volt));
varNoise=zeros(1, length(probe_eater_volt));
for k=1:length(probe_eater_volt)
    bgNoise(k) = mean(result{k}.spectrum.dev3543.demods(4).sample_x_fft_abs_avg{1}.value(bgRange));
    varNoise(k) = sqrt(var(result{k}.spectrum.dev3543.demods(4).sample_x_fft_abs_avg{1}.value(bgRange)));
end
subplot(2, 1, 1)
loglog(probe_eater_volt, probe_power, 'rd-');
xlim([0.01 1]); ylim([0.1 10]);
xlabel('probe eater voltage (V)','Interpreter','latex');
ylabel('probe power (mW)','Interpreter','latex');
grid on

subplot(2, 1, 2)
errorbar(probe_eater_volt, bgNoise, varNoise, 'bo-')
set(gca,'xscale','log');set(gca,'yscale','log');
xlim([0.01 1]); ylim([1e-4 1.2e-3]);
xlabel('probe eater voltage (V)','Interpreter','latex');
ylabel('background noise $\times 100$ (V)','Interpreter','latex');
grid on
%%
save('D:\code\AutoMeasurement\example\DemodPhase\ProbeNoise.mat', 'result');